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Abstract 

Non-Euclidean plates are a subset of the class of elastic bodies having no stress-free 
configuration. Such bodies exhibit residual stress when relaxed from all external 
constraints, and may assume complicated equilibrium shapes even in the absence of 
external forces. In this work we present a mathematical framework for such bodies 
in terms of a covariant theory of linear elasticity, valid for large displacements. We 
propose the concept of non-Euclidean plates to approximate many naturally formed 
thin elastic structures. We derive a thin plate theory, which is a generalization of 
existing linear plate theories, valid for large displacements but small strains, and 
arbitrary intrinsic geometry. We study a particular example of a hemispherical plate. 
We show the occurrence of a spontaneous buckling transition from a stretching 
dominated configuration to bending dominated configurations, under variation of 
the plate thickness. 
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1 Introduction 

Elasticity theory, in its most fundamental formulations, describes the statics 
and dynamics of three-dimensional (3D) elastic bodies. Such "fundamental" 
models are extremely complex, due to both high dimensionality and nonlinear- 
ity. This intrinsic complexity has motivated over the years the development of 
simplified, or reduced models of elasticity. In particular, models of lower spatial 
dimension have been developed to describe the mechanics of slender bodies, 
such as columns, shells and plates. These models are based on various approx- 
imations, such as lateral inextensibility, small deflections and small deforma- 
tions. In particular, the Kirchhoff-Love assumptions [1] allow the derivation 
of reduced two-dimensional (2D) theories of plates. The Foppl-Von Karman 
(FVK) plate equations are one of the successful reduced descriptions of plates 
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mechanics. It expresses the elastic energy of a deformed elastic plate as a sum 
of stretching and bending energies of a 2D surface. The stretching energy, 
which accounts for in-plane deformations, is linear in the plate thickness, h. 
The bending energy, which depends on the curvature of the deformed plate, 
is cubic in h. Other reduced 2D theories usually bear the same structure, i.e., 
their energy is given by the sum of a stretching term and a bending term [2] . 
The validity of the dimensional reduction from 3D to 2D models, based on the 
Kir chhoff- Love assumptions, has been the subject of many scientific disputes 
[3]. Recently, the FVK theory has been derived from a 3D elastic theory by 
means of an asymptotic expansion The stretching and bending terms in 
the FVK theory have also been derived as two different vanishing thickness 
F-limits of the 3D elastic energy [5]. 

2D elastic theories distinguish between two types of thin bodies: plates and 
shells. Plates are elastic bodies that bear no structural variation across their 
thin dimension, and possess a planar rest configuration. Shells are elastic bod- 
ies that bear structural variations across their thin dimension, and as a result, 
possess a non-planar rest configuration. In both cases the postulated existence 
of a stress-free, rest configuration is of paramount importance. 

Recent technological developments have extended the range of mechanical 
structures that can be engineered and constructed. Plates of nanometer scale 
thickness can be manufactured [12], responsive nano-structures are being de- 
veloped [THpHj . and the use of shape memory materials that lead to large 
shape transformations has been extended [15]. In addition, the application 
of mechanics to biological systems, such as in the study of plant mechanics 
and motility [16] and the study of mechanically induced cell differentiation 
[17] . is a rapidly developing field. Such developments have renewed the inter- 
est in elasticity. Several recent theoretical works have focused on the onset 
of various mechanical instabilities and the scaling of the generated patterns 
[T2]IT8] . and other thoroughly analyzed the assumptions underlying some of 
the dimensionally reduced models [5]. 

The modeling of growing elastic bodies is an area in which current theories of 
elasticity face difficulties. Growing tissues, such as leaves, exhibit very com- 
plex configurations even in the absence of external forces [6] . Although leaves 
(and many other growing tissues) are relatively thin (compared to their lat- 
eral dimensions), there are no reduced 2D elastic theories that model their 
shaping mechanisms. Another class of systems for which current theories do 
not apply are elastic bodies undergoing irreversible plastic deformations. The 
main difficulty in applying elasticity theory to growing bodies, or elastic bodies 
having undergone plastic deformations, is their lack of a stress-free configu- 
ration. Specifically, in most models, the elastic energy density of a deformed 
body depends on the local elastic modulus and the strain tensor. The latter 
is defined by the gradient of the mapping between a stress-free configuration 
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and the deformed configuration. It can be shown, for example, that a general 
growth process of an elastic material leads to a body that has no stress-free 
configuration, thus exhibiting residual stress in the absence of external loading 

m- 

To formulate an elastic theory for bodies that do not have stress-free config- 
urations, one needs an alternative definition of the strain tensor. At present, 
certain 3D formulations use the concepts of virtual configuration [Sf^ and in- 
termediate configuration [T9f20] to describe natural growth processes as well 
as plastic deformations leading to residual stress. The growth process in these 
theories is decomposed into a growth step, which maps a stress-free configu- 
ration into a virtual configuration, and an elastic relaxation step, which maps 
the virtual configuration into an elastic equilibrium configuration that contains 
residual stress. These theories use a multiplicative decomposition of the defor- 
mation gradient into an elastic and a plastic part. Other theories decompose 
the strain tensor additively [2T] . 

In the current work, we focus on the elastic response of the body after its "rest 
configuration" has been modified either by growth, or by plastic deformation. 
We do not consider the thermodynamic limitations on plastic deformations 
(which are not relevant to naturally growing tissue). We assume that the dis- 
torted "rest configuration" (or virtual configuration) is a known quantity. If an 
elastic body is capable of assuming the virtual configuration, then there exists 
a stress-free configuration, which is unique; the solution to the elastic problem 
is then trivial. If, however, no elastic body can assume the virtual configura- 
tion, then no stress-free configuration exists, and we face a non-trivial problem 
which exhibits residual stress. We term such bodies as "non-Euclidean" be- 
cause their internal geometry is not immersible in three-dimensional Euclidean 
space. 

We consider now two model examples of elastic structures that belong to 
the class of systems we have termed non-Euclidean plates, and discuss qual- 
itatively some of their properties. Consider an elastic square slab of lateral 
dimensions 2L, and thickness h. Suppose we cut out from it a square segment 
of dimension L, leaving out a U-shaped structure (see Figure [TTj i) . Next, the 
square is replaced by a trapezoid that has three edges of equal length L, and 
a fourth edge of longer size L'. Of course, the trapezoid is too large to fit in 
the square slot. Suppose, however, that we forcefully insert the trapezoid into 
the slot, gluing its three sides of length L to the corresponding edges of the 
U-shape. As a result, the U-shape will slightly open, whereas the trapezoid 
will experience compression. This plane-stress configuration is shown schemat- 
ically in Figure If the plates are sufficiently thin, the trapezoid is unable 
to sustain the compression and buckles out of plane to form a shape qualita- 
tively described in Figure 
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Fig. 1.1. Schematic illustration of an unconstrained plate exhibiting residual stress, 
(a) The two elements composing the plate are shown side by side, (b) As the red 
trapezoid is too large to fit into the square opening, it is compressed, (c) For a plate 
sufficiently thin, the induced compression exceeds the buckling threshold, and the 
trapezoid buckles out of plane . Note that there are many shapes that preserve all 
lengths along the faces of the plate, yet they cannot be planar. 

We note the following points for this toy problem: 

(1) The three dimensional metric that describes the rest lengths of the com- 
pound body (U-shape plus trapezoid) is continuous. 

(2) If denotes the vertical coordinate (say, the distance from the bottom 
face), then all = const surfaces are identical. It is this property that 
causes the body to remain flat (for sufficiently thick samples), and will 
later be used to rigorously define non-Euclidean plates. 

(3) The body exhibits residual stress in the absence of external constraints: 
in Figure |l.lf ) the body is in a state of non-trivial plane-stress, identical 
for all x^ = const sections. In the buckled state (Figure [lT]:) symmetry is 
broken. The upper surface is longer then the lower surface, hence at least 
one of them must be strained. It may easily be shown that the compound 
body has no unstressed configuration. 

(4) The problem is purely geometric: As both pieces (the confining U and 
the trapezoid) are made of the same material, the stiffness of the material 
(Young's modulus) has no effect on the equilibrium shape, and we expect 
to see the same behavior for metals and rubbers (as long as the strains 
are sufficiently small and the stresses are below the yield stress). 

(5) The toy problem presented here may easily be solved numerically us- 
ing commercial software (In fact, a very similar problem was addressed 
experimentally and analytically in [21] )• The treatment used for solving 
such problems is limited to discrete geometric incompatibilities: two (or 
more) regular elastic problems that are coupled through their bound- 
ary conditions are solved simultaneously. Plastic deformations and non- 
homogeneous growth processes, however, cannot be mapped into such 
discrete geometries. 

Recent experiments in torn plastic sheets [22] and environmentally responsive 
gel discs [23] have attracted attention to a specific class of non-Euclidean elas- 
tic bodies: thin bodies whose shaping mechanism is essentially two-dimensional. 
Growing leaves display such behavior, as their growth is believed to be nearly 
homogeneous across their thin dimension, and inhomogeneous in the lateral di- 
mensions. The gel discs reported in [23] mimic a growing thin 3D body shaped 
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by a 2D growth process. In these experiments initially flat stress-free objects 
shrink according to a pre-determined chemical gradient in their composition. 
The shrinking is homogeneous across the thickness, but inhomogeneous in the 



lateral directions (see Figure 1.2 for an example). The resulting body shows 
no structural variation across its thin dimension, yet the lateral equilibrium 
distances, specified by the differential shrinking, define a 2D non-Euclidean 
metric tensor. Thus, they cannot be preserved in any flat conflguration of the 
disc. Such bodies may not be considered as plates (due to their non-planar 
intrinsic geometry), nor as shells (as there are no structural variations across 
the thin dimension). We name such bodies non- Euclidean plates. 

The conflgurations of non-Euclidean plates in the absence of external forces 



are not flat (Figure 1.2: and 1.2 i), and may exhibit multi-scale, and fractal- 
like conflgurations [22 23j. Finite element simulations devised to describe such 
bodies pUIITT] . were able to obtain such multi-scale configurations as energy 
minima. In both computational and theoretical works, it was assumed that 
the elastic energy can be written as a sum of bending and stretching terms. 
The bending was measured with respect to a locally fiat configuration (as in 
the FVK plate model), and the stretching was evaluated with respect to a 
reference 2D metric tensor. None of these works, however, was backed up with 
a theoretical justification for such assumptions. 




Fig. 1.2. An initially flat disc shrinking differentially, (a) The peripheral areas (light 
grey) shrink significantly, while the center of the disc (dark grey) shrinks moder- 
ately, (b) In order to accommodate the center of the disc within the "too short" 
peripheral ring, the plate must buckle out of plane. While the newly prescribed 
lateral lengths are satisfied on average (over the thickness), the symmetry breaking 
causes the upper surface to be tensed, while the lower surface is compressed, (c-d) 
Experimental realization of non-Euclidean plates, using environmentally responsive 
gels as described in [23]. The differential shrinkage prescribes a two dimensional 
geometry of constant positive Gaussian curvature K = O.llcm"^. The thicknesses 
of the gels are he = 0.75mm and hd = 0.6mm. 

In the present work we derive a reduced 2D elastic theory for non-Euclidean 
plates and discuss their characteristics. The derivation starts from a model of 
a 3D covariant "incompatible" elasticity, that is, a model for 3D bodies whose 
intrinsic metric cannot be immersed in a 3D Euclidean space. We advocate 
that the common definition of strains with respect to a stress-free configu- 
ration is too restrictive. Instead, strains can be measured with respect to a 
reference metric tensor, which is not necessarily immersible in 3D Euclidean 
space (incompatibility). When the strain tensor is defined with respect to a 
metric tensor, growth (or any other metric prescription) is naturally decoupled 
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from the elastic relaxation. The second Cauchy-Piola stress tensor (which is 
linear in the strain for small strains), may be written explicitly in terms of the 
difference between two metric tensors. In such a formulation residual stress 
appears inevitably as a result of the lack of immersibility. . 

We apply this formulation to thin elastic plates, using the Kirchhoff-Love as- 
sumptions. When applied to ordinary plates, our theory coincides with the 
Koiter plate theory [2] . As in the FVK and Koiter theories, the energy of the 
plate is a sum of stretching and bending terms. The bending term is cubic in h 
and quadratic in surface curvatures. The stretching term is linear in h and de- 
pends on the difference between the 2D metric tensor of the configuration and 
the reference metric (in [in] it was termed "target metric") The covariant 
elasticity formulation, together with the bending term measures deviations 
from a flat configuration, while the stretching term measures deviations from 
the 2D reference metric (which may be non-fiat). The resulting model is simple 
to use, and has an intuitive structure, which clarifies the underlying physics. 
We end this paper with an application of the theory to a simple case of a 
hemispherical plate. 



2 Theoretical framework: covariant linear elasticity theory 

In this section we derive the energy functional of a three-dimensional elastic 
body as a function of its metric using general curvilinear coordinates. We will 
show that the energy functional takes the following form, 

E{g) = J w{g)\J\g\dx^dx^dx^ w = ^A'^'^'-EijEki, 

where we use the Einstein summation convention and 

A^^ki = Xg^Jf^ + ^ {^g^>^-g^i + fgJ^) e,, = ^ (g,, - g,,) . (2.1) 

Here gij is the metric tensor, cjij is a symmetric positive-definite tensor, which 
we term the reference metric, and X,fi are elasticity (Lame) constants; for 
tensors | ■ | denotes the determinant. This energy functional neglects terms 
that are of order higher than quadratic in Sij, which is the deviation of the 
metric from the reference metric. For bodies which possess a stress-free config- 
uration, g may be called the rest metric and must comply with six additional 
differential constraints (the vanishing of the Ricci curvature tensor). Precise 
definitions will be provided in the following subsections. For a thorough treat- 
ment of bodies that have a stress-free configuration, the reader is referred to 
the recent introductory book by Ciarlet p5], which contains the mathemat- 
ical background to the subject. A similar treatment, which we consider as a 
starting point for our generalization, can be found in Koiter |2] . We derive the 
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energy functional in a slightly different manner, yet we try as far as possible 
to use the notations of [2J, later adopted in [2S]. 

2.1 "Incompatible" covariant three-dimensional elasticity 

When a body (a compact domain Q C M'^) is endowed with a regular set of 
material curvilinear coordinates x = it is also endowed with an 

induced metric tensor. Specifically, if r denotes the mapping from the domain 
of parametrization, V C M.^, into Q (we call r the configuration of the body), 
then the endowed metric is Qij = diV-djV. Here and below we use roman lower- 
case letters ,i,j,... for indices {1,2,3}; the operator di denotes the partial 
derivative with respect to x*. Any deformation of the body (carrying the co- 
ordinates along with every material point) will result in a different metric 
tensor. A rigidity theorem states that if the induced metrics of two configu- 
rations r{x) G Q and r{x) G Q satisfy gij{x) = gij{x) for every x eV, then 
the two configurations can only differ by a rigid motion (a uniform translation 
and a rigid rotation). Thus, the metric (provided that it is immersible in M.^) 
uniquely defines the physical configuration of a three-dimensional body. 

Our main postulate, which may be viewed as a modification of the hyper- 
elasticity principle originally formulated by Truesdell [26j, is: 

The elastic energy stored within a deformed elastic body can be written as 
a volume integral of a local elastic energy density, which depends only on 
(i) the local value of the metric tensor, and (ii) local metrial properties 
that are independent of the configuration. 

The tensors that characterize the material and the body — the elastic ten- 
sors — contain all the information about the elastic moduli and the intrinsic 
geometry of the body. Truesdell's hyper-elasticity principle is formulated in 
terms of the strain tensor, which requires the existence of a stress-free refer- 
ence configuration. In contrast, our postulate is formulated in terms of the 
metric tensor. This obviates the need of a rest configuration, hence allows for 
residual stress. 

Let w be the energy density per unit volume. The total elastic energy is 



Our postulate states that the function w depends on the metric g and on the 
coordinates x (through the elastic tensors), i.e. w = w{g,x). We make the 
following additional assumptions: 




(1) wig,x)>0. 
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(2) For every x E V there exists a unique metric g = g{x) such that 
w{g{x),x) = 0. We call g the reference metric. 

In the present work we consider the reference metric (7 to be a known quan- 
tity, whereas the unknown is g, the "actual" metric of the configuration. It 
turns out to be more convenient to define the energy density per unit volume 
with respect to the volume element induced by the reference metric. We there- 



fore define w = y\g\/\g\w as the new energy density. Note that the previous 
assumptions on w carry over to w, i.e. 

w{g,x)>0, w{g,x) = ^ g = g. 

If we additionally assume that w{g, x) is twice-differentiable with respect to g 
in the vicinity of g, then for small deviations of the metric g from the reference 
metric g our assumptions imply that 

w = lA'''''e,jeki + 0{6^), 

where 

^ij = 2 ~ 9ij ) 

is the deviation of the metric from the reference metric, and A*-'^' can depend 
on g but not on g. 

Note that if there exists a rest configuration {g is an immersible metric), then 
we may choose the coordinates x to be the standard Cartesian coordinates 
on the undeformed configuration, thus setting (jij = 6ij. In such case we may 
define the displacement vector u = r — x to obtain 

e = ^{g-I) = ^{{VrfVr - I) = \{Vu + (Vu)^ + {Vu)^Vu), 

where {'Vr)ij = dvi/dxj. We therefore identify e as the Green-St. Venant 
strain tensor. The Frechet derivative of the energy density w with respect to 
e is the contravariant second Piola-Kirchhoff stress tensor |25] 



For small strains we only need to determine the rank-four contravariant elas- 
ticity tensor A'^^^K Regardless of what g is at any given point p G fi, we may 
always choose a re-parametrization x' such that the reference metric with re- 
spect to the new (local) system of coordinates satisfies g'^j = 6ij at p. If the 
medium is isotropic, then the tensor (^A'y^^'' at p is isotropic in the Cartesian 
coordinates x', hence must be of the form 

^^>yjkl ^ ^^ij^kl ^ ^^^ik^jl ^ ^il^jk^ ^2.3) 
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for some constants A and /i [25] . For a body with a reference rest configuration, 
we may identify these constants as the Lame coefficients. 



It remains to transform the contravariant tensor A', defined on the local Eu- 
clidean coordinates x', back to the original curvilinear coordinates x using the 
transformation rules for tensors, 

j^nmpg ^ ( )« (A-1)™(A-1)P (A"! ^ (2.4) 

where A = dx' / dx is the Jacobian of the transformation (see Appendix |A]) . 
As the strain tensor transforms with the jacobian 

gij - g,j = 2ei, = 2K'lk%i = A'^^g'^ - Sm) = g^j - Af A;.^,, 

we obtain that AfA^-5fc; = (jij. Since all the orientation-preserving Cartesian 
coordinate transformations differ only by a proper orthogonal rotation, this 
equation holds independently of the particular local Cartesian set x'. The 
only implication of this calculation is that g must be symmetric and positive- 
definite, i.e. it is indeed a metric. Yet, this metric is not required to be im- 
mersible in M^, which is why we refer to our theory as "incompatible" elasticity. 

If we now define the reciprocal reference metric by g^^gui = ^1, and sub- 



stitute (Q in using the fact that {A'^Y'^iA-y'dki = we obtain 

expression (2.1) for the energy density. As described in Appendix |A| differen- 
tiation and the lowering and raising of indices are both defined with respect 
to the reference metric. It should be emphasized that Al and 5ij are not ten- 
sors in the sense defined in Appendix [A] {5ij is Kronecker's delta and not the 
lowered- index unit tensor). Moreover, given a metric g^j there exists a recip- 
rocal metric tensor {g~^y^ which is a contravariant tensor of rank two and 
satisfies {g~^y-^gjk = ^1, however it is not obtained by raising the indices of 
gij, i.e. {g~^y-^ 7^ g^^Q^^Qki = g^'' ■ The reference metric is the only tensor for 
which the inverse is obtained by raising both indices. 

The equations of elastic equilibrium are obtained from the energy functional 
by a variational principle. We express the energy as a functional of the metric 
tensor, g, yet variations of g must take into account that its components satisfy 
six differential constraints, which are the vanishing of the Ricci curvature 
tensor. Alternatively, we may vary the configuration r, in which case the 
induced variation in g trivially satisfies the six constraints. Thus, 

5E = [ ^^5eij\f\g\dx = ^ [ S''^ 6 gij \f\g\ dx 
Jv dSij * Jv * 

= j S'^^diV ■ dj5r ^J\g\dx. 

Integrating by parts, and using the fact that 

djdkV = r^i^diV, 
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where 



jk 



lid {djgki + dkgji - digjk 



are the Christoffel symbols associated with the configuration r, we obtain after 
straightforward algebra the following boundary value problem, 



v,s'^ + (r 



I 

jk 



T),)S^' = 
: on dV, 



in V 



(2.5) 



where 



pi 



\g'\djgki + dkgji - dig^t] 



are the Christoffel symbols associated with the reference metric, Uj is the unit 



normal (in 



to dV, and 



M^|\^\S'') + ^kS''' 



is the covariant derivative with respect to the reference metric (see Appendix 
A|. As the elastic body is immersed in the six independent components of 
the symmetric Ricci curvature tensor of the metric g 



Rii = ^{9 {dkdi9ij - dkdjgii + djdigki - didigkj] 



kj'- li 



(2.6) 



must all vanish. The three equations (2.5) together with the six immersibility 



conditions for g (2.6), form a set of nine equations, for the six unknowns in g. 



There are two possible ways to resolve this seemingly over-determination. The 
first is by noticing that the six independent components of the Ricci curvature 
tensor satisfy differential relations: their derivatives are related through the 
second Bianchi identity. The second way of resolving this issue is by identi- 
fying the immersion r as the three unknown functions, in which case the six 



equations in (2.6) are solvability conditions for the PDE (2.5). However, as 



the equations in r are of higher order we need to supply additional conditions, 
namely set the position and the orientation of the body, in order to obtain a 
unique solution for r. 



Eq. (2.5) is our fundamental model for three-dimensional elasticity. The only 



(yet fundamental) difference with standard models of finite displacement elas- 
ticity is that the reference metric does not necessarily have an immersion in 
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3 The elastic theory of non-Euclidean plates 



We define a plate as an elastic medium for which there exists a curvilinear set 
of coordinates in which the reference metric takes the form 



gu 9i2 



9ij = 



921 922 

1 



where d^gij — 0. 



(3.1) 



A plate is called even if the domain P C of the curvihnear coordinates can 
be decomposed into P = 5 x [— |, |], where 5 C and h is constant. Thus 
an even plate is fully characterized by the metric of its mid-surface — 0. 
Let 



dA 



fi'llfi'22 



denote an area element on the mid-surface, and A = dA be the total area of 
the mid-surface. An even plate will be called thin if /i ^ \/A. A plate will be 
called non-Euclidean if the Ricci curvature tensor of its reference metric does 
not vanish. An equivalent condition is that the mid-surface (considered as a 
two-dimensional manifold) has a non vanishing Gaussian curvature. A non- 
Euclidean plate has no immersion with zero strain in M^, i.e. the equilibrium 
state of a non-Euclidean plate must be a frustrated state exhibiting residual 
stress. This statement is rather intuitive: If the plate fully complies with its 
given two-dimensional metric, then it must assume a three-dimensional form 
that violates the invariance along the thin direction. If, on the other hand, it 
remains planar, then it cannot comply with a non-vanishing Gaussian curva- 
ture, hence it must contain in-plane deformations. 



3.1 The reduced energy density 



Although thin plates are three-dimensional bodies, one would hke to take ad- 
vantage of their large aspect ratio and model them as two-dimensional surfaces, 
thus reducing the dimensionality of the problem. Ideally, one would hope to 
obtain a reduced two-dimensional theory as an assumption-free small-/i limit 
of the three-dimensional theory. Unfortunately, such an analysis is still lacking, 
and one must introduce additional assumptions. We adopt the Kirchhoff-Love 
assumptions regarding the structure of the configuration metric g. The stan- 
dard formulation of the Kirchhoff-Love assumptions is: 

(1) The body is in a state of plane-stress (the stress is parallel to the deformed 
mid-surface) . 
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(2) Points which are located in the undeformed configuration on the normal 
to the mid-surface at a point p, remain in the deformed state on the nor- 
mal to the mid-surface at p, and their distance to p remains unchanged. 

The first assumption may be reformulated as 

S'^ = 0. 

In our case, where no reference configuration exists, the second assumption 
may be rewritten as 



9ij 



1 



or equivalently e^s = 0, 



where following [25)12] Greek indices ... assume the values {1,2}. It is 
important to note that the assumptions S''^ = and = represent two 
different elastic problems — plane-stress versus plane-strain problems respec- 
tively. The two stand in contradiction for all A 7^ 0. As a result, the two 
assumptions do not "commute", i.e. the order in which the two assumptions 
are applied is crucial. The key assumption is the first one, S^^ = 0. It states 
that most of the elastic energy is stored in lateral (in-plane) deformations of 
the various constant-x^ planes. Estimates of deviations from this assumption 
may be found in [27]. Let ki and /c2 be the principal curvatures of the mid- 
surface, kmax = max(fci,A;2), and let L be the smallest lateral length scale 
appearing in the elastic equilibrium. It may be shown that the plane-stress 
approximation holds for 

kmaxh -C 1, and h <^ L. 

The second assumption, e^a = 0, is introduced only after we already have a 
reduced energy density, containing only plane-stress contributions. It deter- 
mines the actual three-dimensional configuration the body assumes and the 
variation of the plane-stress along the thin dimension. It enables us to relate 
the elastic energy density to geometric properties of the midplane which is 
considered as a two-dimensional surface. Following [2] we denote by 7 the 
maximal plane-stress of the midplane and note that adding terms of orders 
7^, hkmaxl and h^k"^^^ to the energy density would not modify the order of 
the approximation. Thus the second assumption may be considered as a sub- 
sidiary assumption, used to bring the elastic energy density to the simplest 
consistent form. Although the assumptions are physically plausible, reducing 
the three-dimensional energy functional into a two-dimensional functional by 
means of F-convergence would set the current theory of firmer grounds. 

We now exploit the modified Kirchhoff-Love assumptions to derive a reduced 



two-dimensional model. Combining (2.2) and (2.1) and using the tensorial 
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rules for raising indices we get 

From the first assumption, 5*^^ = 0, and the fact that — ~^ ^3 
^33 = ^3 = follows that 



We use (3.2) to rewrite the energy density (2.1) only in terms of the two- 



dimensional strain, 



or equivalently 

w = lA'^P^'e^pe.s, = 2/i {j^^r^d^' + ^'f'^ ■ 

Note that as we contract the tensors A and A with symmetric tensors we 
only retain their symmetric part. So far we have only used the first of the 
Kir chhoff- Love assumptions. 

We now use the second assumption to express the energy functional as a 
two-dimensional integral over the mid-surface, by integrating w over the thin 
coordinate x^. As 5^33 = d^r ■ d^r = 1 and gas = daV • d^r = 0, we identify 
d^r = N as the unit vector normal to the constant-x^ surfaces. Moreover, 
it can be shown that d^d^r = 0, implying that N = /^/(x^,^^) is the unit 
normal to the mid-surface, and dsd^d^ga/s = 0. 

The most general form of the metric is therefore given by 

gap = dapix^, X^) - 2X^ X^) + (x^)^Ca/3(x\ X^) . (3.3) 

The tensors a, b, c can be identified as follows: we define the mid-surface 

R{x\x^) = r{x\x^,0), 

and note that 

d3gap\x3=o = [dsdaV ■ OpT + d^T ■ d^dpr]^^^^ = -2dadpR ■ N 

and 

dsdsga^l^s^Q = 2da,N ■ dpN, 
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which shows that a,b,c are the first, second and third fundamental forms of 
the mid-surface, i.e. 

aai3 = daR ■ df^R bafs = dadf^R- N = d^N ■ dpN = {a~^y^baybfS5. 

(3.4) 



A metric of the form (|3.3|) with a,b,c given by (|3.4|) corresponds to a three- 

(3.5) 



dimensional configuration of the form 

r(a;\ x^ x^) = i?(x\ x^) + x^iV(a;\ 



Having deduced the dependence of the metric in (3.3), we may integrate the 
energy density over the thin dimension. 



U)2D 



2 I h 

2 



which reduces to 



W2D = -A^'^'^'eiy^f + -^"^^'^ {b^pb.s + el''^{a-Tb,,bs.) + 0{h' 



where e 



2-0 _ l/„ 



24 

(ja/s) is the strain evaluated at the mid-surface. Omitting 
terms of order five and higher in the thickness h, and neglecting e with respect 
to the unit tensor yields the final form of the reduced two-dimensional energy 
density, 

W2n = ^A-^''eiy,f + ^^A-''Xpb,s, (3.6) 

where 

1 + u \1 — u 

We have introduced here the physical constants Y (Young's modulus) and z/ 
(the Poisson ratio), defined by 



2/i 



Y 



and 



A 



1 + u 2/i + A 1-u 

The total elastic energy is obtained by integration over the mid-surface 



We identify the two terms in (3.6) as stretching and bending terms, respec- 



tively, and write the total energy as 



E = hEs + h^EB, 



where 



E. 



/ ws \/\g\ dx^dx"^ Eb = / wb ylgldx^dx"^ 
Js ^ Js ^ 
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and 



Comment; 

1. The quantities Es and are called the stretching and bending contents 
(measures for the amount of stretching and bending that do not vanish in the 
limit — *• 0), and ws and wb are their respective densities. By application 
of the Cayley-Hamilton theorem, the density of the bending content can be 
rewritten in the form 

24(1 + 1^) ^' \g\l 

2. A two-dimensional configuration has zero stretching energy if and only if 
(^ap = Qapi i-G., if the two-dimensional metric coincides with the reference 
metric (such a configuration is an isometric immersion of g). In this case 
(a~^)"'^ = (7"^ and we identify the density of the bending content as the 
density of the Willmore functional [28] 

where K and H are the Gaussian and mean curvatures of the mid-surface. 



3. The total energy (3.7) is a functional of the mid-surface immersion i?, 
i.e., E = E{R). It has two terms: the stretching energy, which scale linearly 
with h, and the bending energy, which scales like the third power of h. The 
equilibrium configuration R* is the one that minimizes the energy functional. 
For thin plates, the total energy is dominated by the stretching term, and we 
expect the equilibrium configuration to have a two-dimensional metric very 
close to the reference metric g. For thick plates, it is the bending energy which 
is dominant, and equilibrium is expected to have a minimal amount of bending. 



3.2 The reduced equilibrium equations 



As in the three-dimensional case, we can derive the Euler-Lagrange equilib- 
rium equations that correspond to the reduced energy functional (3.7) in two 
alternative ways. The first uses independent variations of the six components 
of the symmetric tensors a^^ and 6^^, adding three Lagrange multipliers to 
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impose the three Gauss-Mainardi-Peterson-Codazzi (GMPC) equations: 



N ^ 1 

^1 2 ' \~/^ctp a7 ' ^70- ap - po- J J (39) 



The GMPC equations are the necessary and sufficient condition for and ba/s 
to be the first and second fundamental forms of a surface in M^. It is noteworthy 
that the satisfaction of the GMPC equations is a sufficient condition for the 



immersibihty of a metric of the form (3.3) [25]. Again this mathematical result 



is rather intuitive: If the tensors ttap and ba/3 satisfy the GMPC equations, then 
there exists a mid-surface i?(x^,a;^), for which they constitute the ffist two 



fundamental forms. If such a surface exists then the explicit construction (3.5) 
ensures the existence of an immersion in of the three-dimensional body. 

The second and more natural path is to preform variations in the mid-surface 
R, [25p] . Let us define the reduced two-dimensional stress and moment tensors 
by 

_ 9W2D _ .alS-rS 2D t& _ dW2D _ .^p^s 

' ~ del^p ~ ''^ " db,s ~ 12^ 

Consider then a variation R ^ R + SR. To ffist order in 6R we have 

Sbai3 = dadp5R ■ N + d^dpR ■ 5N = d^dpSR ■ N - Vl^N ■ d^6R, 

where from now on the Christoffel symbols P^^^ are defined with respect to 
the two-dimensional surface R (they are the restriction of P*;, to the indices 
{1,2}). The resulting variation in the energy is 

6E= f (s''^5el^p+m''Hbap)^dx^dx^. 
J s 

Integrating by parts gives the following equation, 

= V, (V/^m"^ + (P^^ - V1p)m'l') - s^^^^p - m'^^Cap 
= (s^^ + m^\a-^r%^,) + (P?^ - f ^^) (s'P + m^^^ia-^VX,) (3.10) 
+ {y,m^^ + (P^, - f ^,)m^'^) {a-y-b,,, 

and boundary conditions: 

= naUpm'^^ 

= np [s'" + {a-^'b^am"^") 

= n, (V„m-^ + (Pf, - Ti,)m'^') , 
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where 



da 



al3 



dp{J\g'\M''^) + f%M^'. 



The three equations (3.10) (in the second equation a = 1,2 is a free index), 
supplemented by the three GMPC equation (3.9), form a boundary value 
problem for a^p and hap as well as an integrability condition for R. 



4 Example: A spherical plate annulus 



4-1 Axially symmetric case 



The reduced two-dimensional equilibrium equations (3.10) are highly nonlin- 



ear equations in the six variables s"^, 771°"^ . A tractable set of equations may 
be obtained if, for example, symmetries are imposed. Let us set = = 9 
(polar coordinates) and consider a reference metric of the following form: 




(4.1) 



In this case, the Gaussian curvature of the mid-surface is = — $rr/$7 where 
we now use subscripts to denote differentiation. Recall that the corresponding 



three-dimensional reference metric Qij given by (3.1) can be immersed in 
only if = 0. 

We seek solutions in the form of a body of revolution 

R{r,9) = (0(r) cos 6',0(r) sin 6', ^/'(r)). 

For such configurations the GMPC equations are satisfied trivially. The first 
and second fundamental forms are given by 



^2 + ^2 



V 



and 



1 



Ja/3 



(bib^^ 



lCr)i then, substi- 



If we define ip^ = 4'rC (which implies that ip^r^^r ~ ipr4>rr - 
tuting the fundamental forms into the two-dimensional energy density (3.6) 
we obtain the following expression for the energy, 



E 



4(1 - 1/2) Js 



W2D^dr, 



(4.2) 
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where 



W2D = hws + h^' 



Wb, 



and 



WS = 2U (02(1 + ^2) _ i) (^2/^2 _ 1) ^ (^(^2(1 ^ ^2) _ 1)2 ^ (^2/^2 _ i)2^) 



are the densities of the stretching and bending contents. Note that the intro- 
duction of C yields an energy density that only includes first-derivatives of 4>, 
and (. 

The minimum energy configuration balances the contributions from both stretch- 
ing and bending terms. Upper bounds on the minimum energy can be derived 

by considering the two extreme cases, which contain no stretching and no 
bending, respectively. Consider first stretch-free configurations, ws = 0, which 
occur when the two-dimensional metric aap coincides with the two-dimensional 
reference metric, gap, i.e., when 



Thus, there exists a unique axially symmetric isometric immersion (however, 
infinitely many non-axisymmetric isometric immersions may exist). The den- 
sity of the bending content of this isometry reduces to 



which is the density ww of the Willmore functional. Integration of this density 
provides a first upper bound on the equilibrium energy. 

Consider next bending- free configurations, wb = 0, obtained if and only if 
C = 0. This implies that V'r = 0, i.e., a flat radially symmetric surface. The 
density of the stretching content reduces to 



Note that there arc infinitely many axially symmetric configurations for which 
the bending content vanishes. Finding the configurations that minimizes the 
stretching energy is equivalent to solving the axially symmetric plane-stress 
problem, which can be achieved numerically. 



= $ 



and 




ws = 2u{<j>l - 1)(07$2 - 1) + {<j>l - If + (0V$2 - 1)^ 
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4-2 Numerical results 



As an example, we consider the case where the two-dimensional reference 
metric Qa/B is that of a sphere, $(r) = sinr, and the domain is an annulus, 

r e [rmin,rmax] C (0,7r/2). 
The stretch-free configuration is a punctured spherical cap and its experimen- 



tal realizations are shown in Figure 1.2 



The minimizer of the energy functional (4.2) was computed numerically for 
the parameters u = 0.5, rmin = 0.1 and rmax = l-l- The elastic modulus Y, 
which is immaterial to the equilibrium shape, was set such that the pre-factor 
7rF/4(l — u^) equals one. As expected, for values of h above the buckling 
transition {hs ~ 0.3) the solution is that of a flat plate, whereas for values of 
h under the buckling transition, the plate is close to spherical. 




Fig. 4.1. Energy scaling in positively curved discs. Total elastic energy, stretching 
energy and bending energy versus the thickness h, for non-Euclidean discs with the 
reference metric and lateral dimensions that are described in the text. All three en- 
ergies are divided by h. The three rendered configurations correspond to (from right 
to left): a flat configuration, a weakly buckled configuration, just below the buckling 
threshold, and a fully buckled, almost isometric configuration (color online). 



In Figure 4.1 we plot the stretching energy (red circles), the bending energy 
(blue crosses) and the total energy (black diamonds) versus the plate thick- 
ness h] all three energies were scaled by 1/h. Except for a narrow transition 
region near the buckling threshold, the total energy is dominated by either the 
stretching energy or the bending energy. As one would expect, the bending 
energy drops to zero above the buckling threshold (large thickness). However, 
below the buckling threshold, as — > 0, the stretching energy drops to zero 
much more rapidly than the bending energy. This last observation is in fact 
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surprising, as naively, one would expect equilibrium to be attained when both 
stretching and bending energy are "equally partitioned" [33] . 



X 10 



Energy 



0.295 0.296 0.297 0.298 
h 



0.299 0.3 



Or 

-0.02- 
-0.04- 
-0.06- 
-0.08- 
-0.1 - 
-0.12- 



Spatial profile 




0.2 0.4 0.6 0.8 



Fig. 4.2. The Buckling transition. The energy variation divided by the thickness, 
(E — EP°'^^)/h (left), and spatial profile (right) near the buckling threshold as a 
function of the thickness. Every profile in the right figure corresponds to an energy 
variation bearing the same color in the left figure. All configurations are fiatter than 



the center disc in Figure 4.1 (note the difference in scales of the axis on the right). 



In Figure 4.2 the spatial profile (a cross-section) of the elastic equilibrium 
configuration is shown. The transition from fiat to buckled configurations oc- 
curs continuously, hence the buckled states, close to the buckling threshold, 
are nearly planar. This supports the validity of theories that assume small 
deflections from a plane (such as the FVK model) for predicting the buckling 
threshold. As the thickness is further reduced, the plate approaches the stress- 
free (isometric) conflguration very fast. The assumption of small deflections 
from a plane fails for such conflgurations. 



The minimal bending content, E^, of the stretch- free conflguration, and the 
minimal stretching content, of the zero bending conflguration yield a 



crossover length scale: he = yE'^'/E'^. Linear analysis about a flat surface 
gives another length scale, the buckling threshold thickness hs- We expect the 
scenario depicted in Figure 4.1 to be valid for bodies in which these two length 
scales are relatively close. However, there are reference metrics (speciflcally, 
hyperbolic), for which all isometric immersions are convoluted, i.e. E^ is very 
large. For such bodies one may obtain he <^ hs- When this occurs, the transi- 
tion region may expand. For such bodies the scaling of the elastic equilibrium 
energy with the thickness will be very different from the one appearing in 
Figure |4.1 
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5 Conclusion 



Natural growth of tissue as well as the plastic deformation of solids are exam- 
ples of local shaping mechanisms of elastic bodies. In general, the local nature 
of such growth processes excludes the existence of stress-free configurations. 
This is the main reason why current elastic theories cannot handle properly 
such shaping mechanisms. In this work we derived a reduced 2D model for 
a class of thin plates with residual stresses, which we named "non-Euclidean 
plates". Such plates are uniform across their thin dimension, but their 2D 
geometry is non-Euclidean. Their complicated 3D configurations cannot be 
obtained from existing 2D models of elasticity. Our derivation is based on a 
covariant formulation of 3D linear elasticity. It does not require the existence 
of a reference stress-free configuration, but only a 3D "reference metric" ten- 
sor, which is determined by the growth. We use this formalism together with 
the Kirchhoff-Love assumptions to derive a 2D energy functional. Like pre- 
ceding theories, this functional decouples into bending and stretching terms. 
The bending term scales like the third power of the thickness and depends on 
surface curvature. The stretching term scales linearly with the thickness and 
increases with in-plane strain, which is nothing but the difference between the 
2D metric tensor of a configuration and the 2D reference metric. Our theory 
is valid for large rotations and displacements and arbitrary intrinsic metrics. 



The numerical results presented in Figure AA_ suggest that in the general case 
there is no equipartition between bending and stretching energies. This in turn 
supports the treatment of very thin bodies as inextensible. Not only the equi- 
librium three-dimensional configuration is dominated by the minimization of 
the "small" bending energy term, but the total elastic energy is dominated by 
it too. The estimate of what thickness should be considered as thin involves 
the introduction of a new length scale he, which is smaller than the buckling 
threshold thickness. The square of this new length scale, h"^, is inversely pro- 
portional to the minimum of the Willmore functional for the prescribed 2D 
geometry. This length scale differentiates between two types of surface geome- 
tries. Surfaces which may be isometrically immersed with a moderate bending 
content, for which he is close to the buckling threshold thickness, will follow 
the shaping scenario and energy profile described in Figures and 4^ Sur- 
faces for which all isometric immersions have high bending contents (as is the 
case for some hyperbolic surfaces) may exhibit very different shaping scenarios 
and energetic landscapes. 



The theory can be further elaborated and generalized to describe a wider range 
of growing bodies. We believe, however, that already in its current stage, it 
is a powerful tool for studying the growth of leaves and other natural slender 
bodies. 
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A Tensors, vectors, scaleirs and the coveiriant derivative 

As our treatment of elastic bodies involves the simultaneous use of two dif- 
ferent metrics, we find it important to provide a brief summary of differential 
geometry in the context of the current work. In the following treatment we 
do not consider the most general setting but only three-dimensional manifolds 
immersed in M^. 

Let the immersed manifold C be the current configuration of an elastic 
body. A global parametrization of is a one-to-one map r : I) ^ Q from 
a domain V C M^. Let r' : V ^ Q he a different global parametrization of 
the current configuration. The composition h = r'~^ o r : T> ^ T>' is called 
a coordinate transformation. The coordinate transformation gradient, often 
denoted by A^ = dx'^ /dx\ is simply the Jacobian matrix of the transformation 
h, i.e. A = dh/dx. The inverse transformation gradient is (A~^)^. = dx^/dx'K 

A scalar is a function $ : — > M. Given a parametrization r : T> ^ Vt, a. scalar 
$ induces a function : P ^ R defined by 0(x) = $(r(x)). Given another 
parametrization r' : V ^ Vl with the coordinate transformation h : P ^ P', 
the relation between the induced functions and (j)' is 0'(x') = 0(h~^(x')). 
By a slight abuse of terminology we also call the functions and 0' scalars 

A vector is a function V from the manifold Q, to the local tangent space of 
the manifold which in our case is M^, V : Q M"^. Note that we cannot 
perform vector operations on pairs of vectors defined at two different points in 
fl, as they belong to different tangent spaces (or equivalently different copies 
of R^). Given a parametrization we may construct a basis e^ = dr/dx^ for 
each tangent space. With respect to this basis we may write any vector as 
V = V^Bi. The three functions : P ^ R arc called the contravariant 
components of the vector V. Again by an abuse of terminology the triplet 
V\ i = 1, 2, 3 is called a contravariant vector. It is easy to prove that under a 
coordinate transformation, a contravariant vector transforms with the inverse 
transformation gradient, — (A~^)^- V'^ , where the left-hand side is estimated 
at a point x while the right-hand side is estimated at the corresponding point 
x' = h(x). 

We next define the dual vector space, namely the space of covariant vectors. 
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However, as there are many ways to define an inner product on the tangent 
space, there are just as many ways to define the dual vector space. The most 
natural inner product is the inner product induced from M^. In such a case, 
we define a dual base e-' by the condition e-' • = 5- , where • is the Euclidean 
product in R^. Any vector in the tangent space may now be decomposed 
with respect to this basis, V = l^e*. The triplet Vi is called a covariant 
vector. Under a coordinate transformation covariant vectors transform with 
the transformation gradient = AjV^'. The inner product in the local tangent 
space induces an inner product on the space of contravariant vectors and the 
mapping of contravariant vectors to their covariant duals by 

V V = V'Wei ■ Bj = QijU'V^ = V%ei ■ = V%6i = UjV\ 

where gij = ■ ej is called the Euclidean metric of Q with respect to the 
given coordinate system. The tensor Qij transforms covariantly in both in- 
dices, i.e. gij = A^A^^f^;. We have identified each contravariant vector with 
a (covariant) vector from the dual space Vi = gijV\ which is called a covari- 
ant vector. The contraction of a covariant and a contravariant vector V^Ui 
yields a scalar. We may choose other inner products on the space of con- 
travariant vectors, leading to different definitions of the dual space. Let cjij 
be a positive definite symmetric tensor, which transforms under a coordinate 
transformation by (jij = A^A^-^j^^ (i.e. covariantly in both indices). The opera- 
tion (, ) : X ^ R given by (U, V) = gijWV^ defines an inner product 
on the space of contravariant vectors. For every contravariant vector there 
corresponds a covariant dual given by Vj — gijV\ The tensor g is called the 
covariant metric on fl. 

Given a parameterized manifold r : P — > Q one may easily prove that the 
gradient of a scalar Vi = di(j) = d(f)/dx^ is a covariant vector. However in order 
to differentiate vectors we need to compare vectors that belong to different 
tangent spaces. To do so we use parallel transport of one of the vectors to the 
point where the other vector is defined. To give only a notion of what parallel 
transport is, we say that it will be transporting the vector along a "straight 
line" , keeping a constant angle between the line and the vector. Both concepts, 
angles between a curve and a vector, as well as "straight lines" (geodesies), 
are defined by the covariant metric tensor. Thus, while the differentiation of 
a scalar is independent of the metric, the differentiation of a vector depends 
on the metric. It may be shown that the parallel transport procedure results 
in the following definition of the covariant derivative. 

ViVj = d,v, - rJ^.Vfc, 

where 

f}fc = Ifidjgki + dkgji - digjk). 
One may verify that VjV^ transforms covariantly in both indices under a 
coordinate transformation. The covariant differentiation of a contravariant 
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vector is given by 

ViV^ = diV^ + flvK 

Note that r*-;^ is not covariant or contravariant in any of its components. 
Henceforth, we will use the term tensors to refer to multidimensional arrays 
for which all indices transform covariantly or contravariantly, thus F is not a 
tensor. One may easily verify that the multiplication or contraction of tensors 
results in a tensor. The differentiation of a tensor should be treated as if the 
tensor is an external product of vectors and apply the covariant derivative 
through the Leibnitz product rule. For example in the two-dimensional case 
we have 

VkM.j = dkMi, - TijMii - TlMji. 

In general, when working with explicit parameterizations we need, in order 
to prove that a certain parameter is a tensor (e.g. a scalar or a covariant 
vector), to prescribe it for all possible parameterizations, and show that it 
obeys the correct transformation rules. This is the case for the current met- 
ric gij = diV ■ djV. It is defined for all possible parameterizations and obeys 
the covariant transformation rules. As the reference metric coincides with the 
current metric (for a local stress- free configuration), we have that g is also a 
rank-two covariant tensor. However some quantities are tensorial by definition, 
for example S^^ = dw/dSij, which is the derivative of a scalar with respect to 
a covariant tensor. For such quantities we may determine their value for one 
(convenient) parametrization, and obtain their value for all other parameter- 
izations through the tensorial transformation rule. This is the case for the 



elastic tensor A*-'^', as may be observed in (2.4) 
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